import string
import copy
import scipy
import Tkinter, tkFileDialog
import numpy as np
import pandas as pd
from scipy.spatial.distance import cdist
import matplotlib.pyplot as plt
import os
import sys
import cPickle
import scipy.io as sio
import matplotlib.units
sys.path.append(os.path.abspath("C:\Users\Scherer Lab E\Documents\GitHub\Python_Data_Analysis"))
from common_functions import *
import half_nanoplate_functions as hnf
import matplotlib_helper_functions as mplhf
%matplotlib inline
#%config InlineBackend.figure_format = 'retina'
#import seaborn as sns
os.chdir("C:\Users\Scherer Lab E\Box Sync\Half-Nanoplate\Figures and Data\Data")
import cPickle
f = open("force_distributions_expt.pkl", "r")
force_dist = cPickle.load(f)
f.close()
f = open("log_pd_experiment.pkl", "r")
log_pd_ls, edge_dist_ls, theta_ls, edge_theta_ls, avg_r_ls = cPickle.load(f)
f.close()
store = pd.HDFStore("K:\Pat's_Projects\ParticleTrajectoryData\half_nanoplate_dynamics_processed.h5", mode='r')
keys = store.index['key']
np_pos = [store.get(v).drop_duplicates(['frame', 'track id']).copy() for v in keys[:-1]]
store.close()
pmf_345 = mplhf.load_matlab_to_dict("pmf_barrier_kernelfit_L345.mat")
umb_samp_pmf = pd.read_csv("emus_pmfs.txt", sep=' ')
def calculated_derivative(x, y, point_shift=7):
"""Calculates the derivative of two Series where one is x and the
other is f(x).
:param x: A pandas Series object of a independent variable
:param y: A pandas Series object of a dependent variable such that
y = f(x).
:param point_shift: The number that points will be shifted for the
derivative. At point_shift=1 the instantanous change in x and y are
calculated. At point_shift=3 the change of x and y are caclulated 3
points in the future."""
if (type(x).__module__, type(y).__module__) == (np.__name__, np.__name__):
change_x = x[point_shift:] - x[:-point_shift]
change_y = y[point_shift:] - y[:-point_shift]
x_points = x[:-point_shift] + (change_x[point_shift+1]/2.0)
else:
change_x = x.shift(point_shift) - x
change_y = y.shift(point_shift) - y
x_points = x + (change_x[point_shift+1]/2.0)
derivative = change_y/change_x
return x_points, derivative
results = []
for l in range(6):
results_dict = {}
# Clean up the US data
x_us = umb_samp_pmf.y*10**3 - 937.5
pmf_us = umb_samp_pmf[str(l)] / (4.11*10**-21)
x_1st_der, der_1st = calculated_derivative(x_us, pmf_us, point_shift=7)
x_2nd_der, der_2nd = calculated_derivative(x_1st_der, der_1st, point_shift=7)
x_spline = np.linspace(-600, 600, 2000)
spline_prams = scipy.interpolate.splrep(x_us, pmf_us, s=25)
spline_y = scipy.interpolate.splev(x_spline, spline_prams)
spline_der_1st = scipy.interpolate.splev(x_spline, spline_prams, der=1)
spline_der_2nd = scipy.interpolate.splev(x_spline, spline_prams, der=2)
results_dict['x'] = x_us
results_dict['pmf'] = pmf_us
results_dict['x_1st_der'] = x_1st_der
results_dict['x_2nd_der'] = x_2nd_der
results_dict['1st_der'] = der_1st
results_dict['2nd_der'] = der_2nd
results_dict['x_spline'] = x_spline
results_dict['spline'] = spline_y
results_dict['1st_der_spline'] = spline_der_1st
results_dict['2nd_der_spline'] = spline_der_2nd
results.append(results_dict)
plt.figure(figsize=(4, 6))
plt.subplot(311)
plt.plot(x_us, pmf_us)
plt.plot(x_spline, spline_y)
plt.title('L='+str(l))
plt.xlabel('Distance from Edge (nm)')
plt.ylabel(r'pmf (kbt)')
plt.subplot(312)
plt.plot(x_1st_der, der_1st, 'o')
plt.plot(x_spline, spline_der_1st)
plt.plot((-600, 600), (0,0), '--k')
plt.xlabel('Distance from Edge (nm)')
plt.ylabel(r'$\frac{d}{dx}\mathregular{pmf_{US}}$')
plt.subplot(313)
plt.plot(x_2nd_der, der_2nd, 'o')
plt.plot(x_spline, spline_der_2nd)
plt.plot((-600, 600), (0,0), '--k')
plt.xlabel('Distance from Edge (nm)')
plt.ylabel(r'$\frac{d^2}{dx^2}\mathregular{pmf_{US}}$')
plt.show()
def nearest_index_in_array(val, array):
return np.argmin(np.abs(array - val))
Calculate where the first derivative crosses 0 for L=0-3 to find the minima and maxima that make up the barrier. From this can calculate the $\Delta U$ and $x_b$ from each pmf. The derivative is calculated from the data directly. Zero crossing point in x is found using linear interpolation between the two points that are nearest to crossing 0. The x values are used in a spline fit of the original data to calculate $\Delta U$.
x_bs = []
del_gs = []
for num, l_dict in enumerate(results):
deriv = l_dict['1st_der']
x_deriv = l_dict['x_1st_der']
zero_crossing = deriv[(deriv.shift(1) * deriv) < 0]
if len(zero_crossing) <= 1:
continue
index = zero_crossing[:2].index
x_crossing = []
# Find where it crosses 0 using linear interpolation
for idx in index:
f = scipy.interpolate.interp1d(deriv.loc[idx-1:idx].values, x_deriv.loc[idx-1:idx].values, kind='linear')
x_cross = f(0)
x_crossing.append(x_cross)
original_indx = [nearest_index_in_array(v, l_dict['x_spline']) for v in x_crossing]
x_b = x_crossing[-1] - x_crossing[0]
org_pmf = l_dict['spline']
del_g = org_pmf[original_indx[1]] - org_pmf[original_indx[0]]
print "L = "+str(num)
print "Minimum at "+str(x_crossing[0])+" nm and "+str(org_pmf[original_indx[0]])+" kbT"
print "Maximum at "+str(x_crossing[1])+" nm and "+str(org_pmf[original_indx[1]])+" kbT"
print "x_b = "+str(x_b)
print "del_g = "+str(del_g)+"\n"
x_bs.append(x_b)
del_gs.append(del_g)
The minima and the maxima in the 2nd derivative are found as they should correspond to the minima and maxima in the original function. A maxima in the 2nd derivative corresponds to a minimum in the original function and a minimum in the 2nd derivative corresponds to a maximum in the original function. Some of these second derivative curves have two maxima near where a minimum is expected in the orginal function. The one chosen is the one with the highest x value. No interpolation was done in this analysis as the max/min was just found on an array (within some bounds) to determine the maxima and minima of the 2nd derivative. From the x values of the maxima and minima, $x_b$ and $\Delta U$ were calculated.
x_bs_2nd = []
del_gs_2nd = []
for num, l_dict in enumerate(results):
x = l_dict['x_2nd_der']
y = l_dict['2nd_der']
maxima = np.argmax(y[np.logical_and(x < -200, x > -350)])
if num <= 1:
maxima = np.argmax(y[np.logical_and(x < -200, x > -500)])
minima = np.argmin(y[np.logical_and(x < 0, x > -100)])
indx_min = nearest_index_in_array(x[minima], l_dict['x'])
indx_max = nearest_index_in_array(x[maxima], l_dict['x'])
pmf = l_dict['pmf']
x_b = x[minima] - x[maxima]
del_g = pmf[indx_min] - pmf[indx_max]
print "\nL="+str(num)
print 'x_b = '+str(x[minima] - x[maxima])
print 'del_g = '+str(pmf[indx_min] - pmf[indx_max])
x_bs_2nd.append(x_b)
del_gs_2nd.append(del_g)
Taken from the first derivative analysis.
fig = plt.figure(figsize=(3.35, 2))
Ls = range(len(x_bs))
ax1 = plt.subplot()
plt.plot(Ls, del_gs)
plt.xlabel(r'$l$', usetex=True)
plt.ylabel(r'$\Delta G^\ddagger$ or $\Delta U$ ($k_B T$)', usetex=True)
ax1_2 = ax1.twinx()
ax1_2.plot(Ls, x_bs, 'r')
plt.ylabel(r'$x^\ddagger$ (nm)', usetex=True, fontdict={'color':'red'})
plt.show()
params = {}
params['ls'] = Ls
params['x_bs'] = x_bs
params['del_gs'] = del_gs
'''Save the componenets of the del_g and x_0 graph'''
import cPickle
pickle_file = open("C:\Users\Scherer Lab E\Box Sync\Half-Nanoplate\Figures and Data\Data\\graph_del_g_x_o.pkl", 'w')
cPickle.dump(params, pickle_file)
pickle_file.close()
Taken from the second derivative analysis.
fig = plt.figure(figsize=(3.35, 2))
Ls = range(len(x_bs_2nd))
ax1 = plt.subplot()
plt.plot(Ls, del_gs_2nd)
plt.xlabel(r'$l$', usetex=True)
plt.ylabel(r'$\Delta G^\ddagger$ or $\Delta U$ ($k_B T$)', usetex=True)
ax1_2 = ax1.twinx()
ax1_2.plot(Ls, x_bs_2nd, 'r')
plt.ylabel(r'$x^\ddagger$ (nm)', usetex=True, fontdict={'color':'red'})
plt.show()